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Abstract 

Laser control of ultrafast double proton transfer is investigated for a two- 
dimensional model system describing stepwise and concerted transfer path- 
ways. The pulse design has been done by employing optimal control theory 
in combination with the multiconfiguration time-dependent Hartree wave 
packet propagation. The obtained laser fields correspond to multiple pump- 
dump pulse sequences. Special emphasis is paid to the relative importance 
of stepwise and concerted transfer pathways for the driven wave packet and 
its dependence on the parameters of the model Hamiltonian as well as on 
the propagation time. While stepwise transfer is dominating in all cases con- 
sidered, for high barrier systems concerted transfer proceeding via tunneling 
can make a contribution. 
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1. Introduction 

Laser pulse design has developed into a powerful method for control- 
ling molecular dynamics [H, 0, 0, 0, 0, 0]. In particular laser control of 
proton transfer reactions in the electronic ground state has triggered con- 
siderable theoretical efforts A straightforward method for laser-driven 
proton transfer between two minima on a potential energy surface (PES) 
is the so-called pump-dump approach adapted from the study of isomeriza- 
tion reactions jl, 0, In Ref. 11] pump-dump control was demonstrated 



for a two-dimensional (2D) model mimicking the situation in malonaldehyde 
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derivatives. The first pump-pulse promotes the system from the localized 
reactant state into a delocalized state above the reaction barrier from where 
it is dumped into the product well by a second pulse. The optimization of the 
pulse parameters can be done by hand, but it was found that optimal control 
theory (OCT) [fjl, 13, 14 1 applied to this problem essentially yields a similar 



pump-dump pulse field [15|, [16(. OCT provides a means for controlling the 
pulse intensity via a penalty factor. Increasing the penalty for strong pulses 
in the case of proton transfer reactions resulted in a change of the mechanism 



from above barrier transfer to tunneling through the barrier [15l. Il6l |. In this 



case the resulting driving field resembles a few-cycle pulse and the associated 
dynamics in isolated and dissipative systems was investigated in some detail, 
e.g., in Refs. [13, [3| (for applications of few-cycle pulses to isomerization 
reactions, see e.g. Refs. [l^, Ho])- Besides these efforts along a quantum me- 
chanical description, there have been attempts to trigger transfer in a double 



minimum PES using classical trajectory-based local control theory [21 



Laser control of double proton transfer (DPT) reactions has received less 



attention. Nishikawa et al. [22[ considered control of stepwise DPT in an 



asymmetrically substituted tetraflouro-porphyrin model using the stimulated 
Raman adiabatic passage (STIRAP) method. Assuming a stepwise transfer 
the description has been reduced to two independent one- dimensional poten- 
tials obtained from intrinsic reaction coordinate calculations. It was shown 
that STIRAP can achieve population inversion related to stepwise DPT on 
a time scale of some tens of picoseconds. Shapiro and coworker 23J have 
investigated DPT in the context of DNA radiation damage and repair for a 
dinucleotide model. They developed a two-dimensional potential energy sur- 
face for DPT comprising a linear reaction coordinate for the concerted DPT 
supplemented by an out-of-plane squeezing type vibrational mode. Employ- 
ing coherently controlled adiabatic passage, detection and repair of DPT 



related mutation has been demonstrated. Finally, Thanopulos et al. [24 
have shown the implications of DPT control for single molecule charge trans- 
fer. Considering a thio-functionalized porphyrin derivative attached to four 
gold electrodes it was found from molecular orbital analysis that the differ- 
ent trans-forms essentially provide orthogonal pathways for electron transfer. 
Laser controlled switching was discussed in a stepwise model comprising two 
bond coordinates for the first step and a single linear reaction coordinate for 
the second one. 

The motivation for the present study has been twofold. First, porphyrin- 
derivatives are not only interesting for charge transfer, but are used in the 
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context of excitation energy transfer as well 25(. The rate for excitation 



energy transfer between two chromophores depends on the relative orienta- 



tion between their transition dipole moments [26|. Since this orientation is 
affected by DPT in the electronic ground state one can, at least in principle, 
envision a laser-triggered ultrafast switch built into an energy transfer de- 
vice. Second, previous studies on controlled DPT considered either stepwise 
or concerted pathways. The decision for a pathway which is dominating in 
an actual system is often not that clear-cut (see discussion of the porphycene, 
a structural isomer of porphine in Ref. |27|) and may depend on the experi- 
mental conditions such as temperature [28j . Furthermore, the mechanism of 
DPT is commonly inferred from thermal rate constants or line splittings in 
absorption spectra. Driving a system with a control pulse, however, is likely 
to establish a pronounced nonequilibrium situation. Here, it is not evident 
that the system will follow the same reaction path as in the quasi-equilibrium 
case. 

In the following Section[2]we will present a two-dimensional model Hamil- 
tonian which comprises concerted and stepwise DPT. This allows us to iden- 
tify their relative importance in dependence on the applied control field. 
Laser control is achieved employing OCT in the frame of a multiconfiguration 
time-dependent Hartree (MCTDH) simulation of the wave packet dynamics 



23, [30|. Respective results will be discussed in Section [3j 



2. Theory 



2.1. Model Hamiltonian 

The DPT Hamiltonian will first be given in terms of single proton transfer 
coordinates x% and x 2 which are assumed to describe the linear translocation 
of the particles between donor and acceptor sites. The simplest form is ob- 
tained by combining two quartic potentials coupled via a bilinear interaction 
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This type of potential has been used extensively by Smedarchina and cowork- 
ers |3 ll . |32| . It has two equivalent global minima (called "trans" in the fol- 
lowing) and two equivalent local minimum (called "cis") which correspond to 
the transfer of a single proton. In eq. ([I]), Uq is the barrier for uncorrelated 
single proton transfer, Xq is half the transfer distance, and g is the dimen- 
sionless bilinear coupling strength. Since we are aiming at an analysis in 
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terms of concerted and stepwise transfer, this potential is more conveniently 
expressed in terms of the symmetric and asymmetric transfer coordinates 
x s = (xi + Xq)/2 and x a = {x\ — Xq)/2, respectively. This transformation 
gives for eq. ([I]) 

U sym (x s , x a ) = 2U + ^ \{g - A)x 2 a - (g + A)x 2 s 



x, 





+ j~{x s x a 6x s x a ) . (2) 
x 

For the purpose of laser control of DPT it is useful to consider the case of 
asymmetric molecules in order to allow for a clear identification of initial 
and final states. Two types of asymmetry can be introduced into the model 
Hamiltonian, that is, with respect to the trans and cis states: 

TT i \ ^trans U , a cis U 

^asym I X s , X a ) X s T" X a yO) 

x x 

where atrans and a C j S are dimensionless parameters characterizing the detun- 
ing between the trans and cis states, respectively. 

We will choose parameters such as to highlight in particular cases of high 
and low barriers for the concerted DPT. The reference case will be the high 
barrier scenario with equivalent cis minima shown in Fig. [1] (for parameters 
see figure caption). Note that we do not address any specific system in our 
study, although the energetics of the PES can be considered to be typical for 
DPT molecules. Here, in particular porphycene derivatives offer a tunability 



of the energetics over a wide range (see, e.g., [331] ). 
The total Hamiltonian is given by 

h 2 ( d 2 d 2 \ 

H = _ 2^ ( ft^2 ~*~ 'dx 2 ' j ^ sym ( Xs, Xa ' ~*~ ^asym^' X «) (^) 

where it is assumed that the moving particles are H-atoms and thus the mass 
for the collective coordinates is m = 2mn- 

For the interaction with the laser field we assume that the permanent 
dipole moment depends only linearly on the coordinates, i.e., 

#field(*) = ~{d a Xa + d s X s )E(t) (5) 

where d a / s is the derivative of the dipole moment with respect to x a / s and 
E(t) is the laser field. For the case of symmetric systems like porphycene the 
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dipole moment changes only along the asymmetric coordinate (for simplicity 
we have chosen d a = l.Oe). This will be used in the reference case. In 
order to investigate the principal effect of a dipole gradient along x s due to 
asymmetric substitution we will also consider a variation along x s . 

2.2. Quantum Dynamics 

The two-dimensional time- dependent Schrodinger equation 



d 

ih—V(x a , x s ; t) = (H + tffieid (t))*(ar , x,l t) 



(6) 



has been solved using the MCTDH approach as implemented in the Heidel- 
berg program package 34j. To this end the two-dimensional wave function 
^(x a , x s ; t) is represented on a grid in terms of a harmonic oscillator discrete 
variable representation (64 points within [-2.5:2.5]a ). The actual propaga- 
tion is performed using the variable mean field scheme in combination with 
a 6-th order Adams-Bashforth-Moulton integrator. Using 20 single particle 
functions per coordinate the largest natural orbital populations have been 
typically on the order of 10 -6 . Selected eigenstates of the time-independent 
Hamiltonian, 0j, are obtained by improved relaxation [35| . 

The shape of the laser field is determined using OCT, that is, the following 
functional is maximized at some final time T 36j, |37| 



J(E,T) = (V(T)\6\V(T)) - k dt 



E\t) 



sin 



(Trt/T) 



(7) 



where O is the target operator and k the penalty factor. The iterative opti- 
mization of the field has been performed by using the implementation within 
the MCTDH program package by Brown and coworkers as detailed in Ref. 
38j. Note that their approach doesn't follow the standard iteration proce- 



dure which is not suitable in the context of MCTDH. Instead an extrapolation 
scheme is used following the suggestion of Ref. 
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The pulses will be characterized in terms of their XFROG trace 



'XFROG 



drE(r)G(T - t)e 



(8) 
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where G(t) is a step-like gate function with Gaussian tails 

It will be assumed that initially the system is in the vibrational ground 
state which is localized in the left part of the potential due to the asymmetry 
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term a trans . The target operator will be taken as the projector onto that 
state which is localized in the right well of the potential, i.e. O = \4>i){<f>i\ 
(see Fig. [3]). In order to elucidate the mechanism of DPT, i.e., stepwise vs. 
concerted, we have defined step-like operators dividing the different regions 
of transfer between the two trans-forms as shown in Fig. HJ Here, transfer 
is counted as concerted if the wave packet passes a narrow range in the 
vicinity of the second order saddle point of the PES. The width of that 
range is, of course, arbitrary and we have used the width of the ground 
state distribution. Furthermore, note that this definition of "stepwise doesn't 
imply the existence of a stable intermediate in the sense of traditional kinetic 
considerations. 

The results reported below have been obtained for different numbers of 
OCT iterations starting with a sin 2 -shaped guess field. The convergence of 
the control functional has been monotonic. The iteration number has been 
chosen such that the change in the control yield between two iterations was 
below 10 -4 . The guess field frequency we have chosen to correspond, for the 
high barrier case, to the transition between the ground state and the first 
excited state along the x s coordinate which is slightly red-shifted as compared 
with the strong IR active first excited state along the x a coordinate. For the 
low barrier case the resonance to the x a coordinate has been chosen. Lowering 
the guess frequency, e.g. by 20 cm -1 , resulted in zero yield, an increase by 
the same amount gave a similar convergence behavior of the control yield. 

3. Results and Discussion 

3.1. Dependence on Pulse Duration 

The OCT optimization problem depends on the final time T at which the 
optimization goal shall be reached. The obvious question is to what extent 
will the OCT field and the associated wave packet dynamics depend on the 
pulse duration. Fig. [2] compares OCT pulses for T = 500 (a) and 1500 fs 
(d). Inspecting their XFROG traces in panels (b) and (e) we notice that 
upon increasing T one obtains a pulse train like structure with an increasing 
number of overlapping subpulses. 

The dynamics can be analyzed in terms of populations of eigenstates of 
the field-free Hamiltonian. A level scheme as well as densities for selected 
states are given in Fig. [31 In total there are 24 states below the barrier for 
concerted DPT. Also shown are the localized initial, <po, and target, <pi, states 
which differ by 38 cm -1 for the chosen a trans . The population dynamics for 
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T = 500 fs is presented in Fig. HI Panel (a) shows that the initial state is 
almost completely depopulated during the first ~180 fs whereas the target 
state becomes populated after ~300 fs reaching a final population of 0.74. 
Since the dipole moment changes along the x a direction only, we would expect 
that the laser pulse excites a wave packet along this direction. Indeed the 
initial excitation is to state 06 which has a node along x a . But already the 
next state which is populated, i.e. 0i6, is of mixed character containing 
excitations along both directions. This mixing is, of course, a consequence 
of the anharmonicity of the potential. Both states are excited by the first 
pulse within about 200 fs. The mismatch between the 0o ~~ ► 06 and 06 — > 016 
transitions amounts to a 65 cm -1 redshift which is covered by the pulse 
spectrum. The latter also shows a slight down-chirp. Due to the delocalized 
nature of states 0i6 and 4>n the first pulse actually excites a superposition 
of these two states, what can be seen from the rise of the population of state 
017 around 100 fs in Fig. HI The second pulse which is centered around 300 
fs acts as a dump pulse transferring the populations according to the scheme 
017 ~~ * $7 — * <fii- in between the two main pulses there is a minor feature at 
~220 fs spectrally located around 200 cm -1 . It is related to the population 
of states which are energetically above the barrier for concerted transfer (see, 
Fig. EJ). In fact the two main pulses excite and de-excite state 024 according 
to 0i6 — > 024 and 024 — > 0i7, respectively. The subpulse centered around 200 
cm -1 switches populations according to 2 4 — > 027 — > 024- This effect is most 
likely not relevant for the DPT control. 

The population dynamics for the T = 1500 fs case is shown in Fig. [5j As 
compared with the T = 500 fs case it is more complex. Only the effect of the 
first pulse is comparable, that is, the excitation scheme 0o — > 06 — > 0ie/0i7 
applies. Some features of the dynamics triggered by the other subpulses are: 
(i) The initial state is repopulated around 600 fs and subsequently depop- 
ulated until about 900 fs. (ii) A population/depopulation is also observed 
for the target state around 1000 fs. (iii) Above barrier states are populated 
directly by transitions from below barrier states. Finally, the population of 
the target state at 1500 fs is 0.91, thus exceeding the one for the 500 fs case, 
i.e. the convergence is faster for the longer pulse. 

Since the reference case contains a dipole moment in x a direction only 
one expects that the laser driven wave packet is bound to move toward the 
product well in a stepwise fashion, i.e. passing through the regions Si and 
S2 in Fig. [TJ Indeed this is the major pathway as shown for different pulse 
durations in Figs. and f. We further notice that the rise of the target 
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state population is closely related to the decay of the probability of being in 
the regions Si and S 2 . For example, for T = 1500fs (Fig. (2U-f) the latter has 
decayed by ~ 1200 fs which coincides with a steep rise of P\. In other words, 
it takes about 1200 fs for the major part of the laser-driven wave packet to 
reach the product well via a stepwise mechanism. Similar arguments hold 
for the T = 500 fs case where this transition takes place around 300 fs. 

Comparing the results for the two propagation lengths one notices a dif- 
ference as far as the dynamics in the range of concerted DPT is concerned. 
Inspecting the densities of the different populated states it is clear that the 
contribution to the concerted pathway is mostly due to the pair 0i 6 and 4>vr- 
A superposition of these states is excited already by the first pulse and the 
subsequent dynamics will be that of a two-level system. The energy mis- 
match between these two states of 51 cm -1 yields a time scale of ~650 fs 
for a round trip between the two minima, i.e. the transfer from the reactant 
to the product well takes about 325 fs. Given the approximate nature of 
this estimate, e.g. due to preparation process, one finds this time scale for 
the motion through the concerted region in the T = 500 fs case. For the 
T = 1500 fs case the time that the wave packet spends in the C region is 
about 900 fs and from Fig. [2]f one notices that there are three maxima of the 
C region population during this period. Since 0i6 and (f>n are energetically 
below the barrier, one can conclude that the contribution to concerted DPT 
is mostly due to tunneling. 

3.2. Dependence on Model Parameters 

Lowering the reaction barriers has a substantial influence on the OCT field 
and the wave packet dynamics. Exemplarily we show the case of T = 1500 
fs in Fig. O First, we notice that, in contrast to the high barrier case, the 
OCT pulse consists of only three subpulses. These pulses drive the wave 
packet essentially via the stepwise pathways as seen in Fig. [6b. This is 
somehow counter-intuitive since one would expect that upon lowering the 
barrier for concerted DPT, tunneling becomes even more effective. However, 
for the present parameters there is essentially only one state below the barrier 
having locally an excited x a character, namely 0s; see Fig. [7J This state is not 
appreciably mixed with other states, say of x s excitation character. Instead 
it is delocalized encompassing both trans as well as the cis regions. The 
population dynamics in Fig. [8] shows that the dominant pathway is indeed 
0o — > 05 ~ *• 0i triggered in a pump-dump like fashion. The spectrally broad 
pulses also excite state 07 which in turn causes an excitation of an above 
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barrier state, <pi2- The latter state clearly will lead to a contribution of the 
concerted pathway, but this time not by virtue of tunneling. However, it's 
effect is rather minor as can be traced from the population of the concerted 
region in Fig. [6b. Finally, we note that due to the sparser energy level 
structure the OCT iteration converges faster and a 99% population of the 
target state is achieved. 

Next we explore the effect of adding an asymmetry to the two cis con- 
figurations for the high barrier reference case. We have chosen a C i S = 0.01 
which translates into an energetic difference of 77 cm -1 between the two cis 
minima. In the symmetric case the populations of the different cis regions Si 
and S2 reflect the oscillation of the wave packet excited in the more stable 
trans well as seen in Fig. [Hh- That means, if the major part of the wave 
packet moves via Si there is a minimum in the population of S2 and vice 
versa. Overall there is no net preference for a pathway. If the S2 region is 
energetically more favorable OCT prefers this way as seen from Fig. [9b. 

Finally, we address the case where there is a nonvanishing dipole moment 
along both directions, exemplarily we have chosen d a /d s = 4. Here, direct 
excitation of the x s coordinate becomes possible and one would expect the 
importance of the concerted pathway to increase. Indeed this is the case and 
for the optimized field the integrated populations for both pathways during 
the time interval up to 1500 fs are obtained to give a ratio of C to t/Stot = 0.144. 
This value exceeds the 0.130 obtained from Fig. [2J. However, the extra gain 
in concerted DPT is small. In fact the OCT pulse populates a state of 
x s excitation character (0 4 ) which is ~36 cm -1 below (pQ. The population, 
however, doesn't exceed 10% and is subsequently also promoted to 0i6- Thus 
the smallness of the gain in concerted DPT is due to the fact that <piQ, which 
plays the dominant role for the concerted DPT, can be reached from the x a 
and x s fundament excitation state 4>q and 04, respectively. That is, there is 
no substantial net effect due to the (competing) excitation of the symmetric 
coordinate. 

In summary, optimal control theory has been used to devise laser fields 
which trigger ultrafast double proton transfer in a simple two-dimensional 
potential supporting stepwise and concerted mechanisms. The optimized 
pulses are of (multiple) pump-dump character, stepwise transfer was found 
to be dominant for all cases considered. Concerted transfer has been shown 
to become relevant if almost degenerate pairs of states exist which are mostly 
localized in the reactant or product well but can be excited simultaneously. 
In order to excite such a superposition state for the situation where the dipole 
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moment changes along the asymmetric coordinate only, the states have to be 
of mixed character with respect to symmetric and asymmetric excitations. 
This situation is likely to be met for high barrier cases which support many 
locally excited states, that attain mixed character with increasing energy due 
to the anharmonicity of the potential energy surface. In such a situation a 
dipole change along the symmetric coordinate does not provide an additional 
means for concerted DPT. For low reaction barriers already fundamental 
excitations of the asymmetric coordinate can be delocalized such as to provide 
a direct pump-dump pathway involving a mostly single intermediate state. 

The present use of a rather simple model Hamiltonian facilitates applica- 
tion to specific molecular systems by virtue of its straightforward parametriza- 
tion. An extension of the latter might be necessary, however, since it is known 
that the effect of heavy atom motions on the proton transfer and the mech- 
anism of DPT can be substantial. In this respect the used OCT-MCTDH 
approach is ideally suited for studying driven multidimensional quantum dy- 
namics. 
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Figure 1: (color online) Two-dimensional PES for the high barrier case (Uq = 2000 cm -1 , 
xq = lao, g = 0.2, atrans = 0.005, and a C i S = 0.0). The energy is given in units of the 
barrier height for concerted transfer (4428 cm -1 ). The energetic difference between the 
two trans minima is 41 cm . The two dashed lines mark the region which is passed 
during stepwise DPT. This doesn't include the shaded area whose width in x a direction 
corresponds to the width of the vibrational ground state and which is assigned to concerted 
DPT. Upper and lower stepwise DPT is marked by Si and S 2 , respectively. In the low 
barrier case (not shown) Uq = 800 cm -1 , giving a barrier height of 1773 cm -1 . The 
barriers for stepwise DPT are 2430 cm" 1 and 975 cm -1 and the energies of the cis minima 
are 821 cm -1 and 330 cm -1 for the high and low barrier case, respectively. 
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Figure 2: (color online) Laser pulses obtained from OCT using (a) T — 500 fs (guess field 
E = 7 mE h /ea ) after 291 iterations, (d) T = 1.5 ps, (guess field E = 0.95 mE h /ea ) 
after 165 iterations (k — 1.5 a.u., uj = 1194 cm -1 ). Panels (b) and (e) show the respective 
XFROG traces for a gate function of width 1000 a.u. and Gaussian tails of width 1000 
a.u. (contours (b) 2.6 to 8.9 in steps of 0.3 a.u., (e) 1.4 to 4.6 in steps of 0.2 a.u. ). Panels 
(c) and (f) give the probability of being in the C and S regions of the PES. 
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Figure 3: (color online) Level scheme for the high barrier case and densities of selected 
eigenstates which get appreciably populated during the laser driven transfer. The lengths 
of the vertical solid bars indicate the range of derealization of the state with respect to 
the symmetric coordinate x s (e.g. 4>o/4>i are localized in the left/right trans minimum, 
02/03 (^2000 cm -1 ) around x s = 0, i.e. in the cis minimum, etc.). The dashed lines 
correspond to the energy of the barriers for concerted and stepwise transfer. The density 
plots cover the range [—1.5 : 1.5]ao along the vertical x a and the horizontal The 
vertical arrows indicate the major pathway for laser driven transfer, the horizontal one 
shows the pathway for concerted transfer. (Note that around 4200 cm -1 there are several 
close lying states.) 
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Figure 4: (color online) Population dynamics for the pulse of Fig. [2^ and the states shown 
in Fig. [3] 
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Figure 5: (color online) Population dynamics for the pulse of Fig. [2ji and the states shown 
in Fig. [31 ( In panel (b) J^and Pn have been offset by 0.6 for clarity.) 
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Figure 6: OCT simulation for the low barrier case, (a) Optimized field for T = 1500 fs 
(guess field Eq = 0.95 mEh/eao, u> — 637 cm -1 , k — 1.5 a.u.) after 60 iterations, (b) 
XFROG traces for a gate function of width 1000 a.u. and Gaussian tails of width 1000 
(contours from 0.5 to 1.7 in steps of 0.1), (c) probability of being in the two cis (S) or in 
the concerted (C) transfer region. 
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Figure 7: (color online) Level scheme for the low barrier case and densities of selected 
eigenstates which get appreciably populated during the laser driven transfer. The density 
plots cover the range [—1.5 : 1.5] ao along the vertical x a and the horizontal (for 
meaning of arrows and vertical bars, see Fig. [3|) 
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Figure 8: (color online) Population dynamics for the pulse of Fig. [6^ and the states shown 
in Fig. El 
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